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Adaptive 2D IGA boundary element methods 


Michael Feischl, Gregor Gantner*, Alexander Haberl, Dirk Praetorius 


Abstract 

We derive and discuss a posteriori error estimators for Galerkin and collocation IGA boundary element methods for 
weakly-singular integral equations of the first-kind in 2D. While recent own work considered the Faermann residual error 
estimator for Galerkin IGA boundary element methods, the present work focuses more on collocation and weighted- 
residual error estimators, which provide reliable upper bounds for the energy error. Our analysis allows piecewise 
smooth parametrizations of the boundary, local mesh-refinement, and related standard piecewise polynomials as well 
as NURBS. We formulate an adaptive algorithm which steers the local mesh-refinement and the multiplicity of the 
knots. Numerical experiments show that the proposed adaptive strategy leads to optimal convergence, and related IGA 
boundary element methods are superior to standard boundary element methods with piecewise polynomials. 

Keywords: isogeometric analysis, boundary element method, collocation, a posteriori error estimate, adaptive 
mesh-refinement 


1. Introduction 


1.1. Isogeometric analysis 

The central idea of isogeometric analysis (IGA) is to 
use the same ansatz functions for the discretization of the 
partial differential equation at hand as for the representa¬ 
tion of the problem geometry. Usually, D is represented 
in computer aided design (GAD) by means of NURBS, hi- 
era rchical s plines, or T-splines. This concept, invented 
HGB05l | for finite element me thods (IGAFEM) has 


proved very fruitful in applica tions HGBOSl |SBTRI2| ; see 
also the monograph CHB09l |. Since GAD directly pro¬ 
vides a parametrization of the boundary dO, this makes 
the boundary element method (BEM) the most attrac¬ 
tive numerical scheme, if applicable (i.e., provided that the 
fundamental solution of the differential operator is explic¬ 
itly known). Is ogeometri c BEM (IGAB EM) has f irst been 
considered in PGK~*~09l | for 2D resp. ISSE+UI^ for 3D. 
While standard BEM with piecewise polyno mials i s well- 
studied in the literature, cf. the monographs SSIlL Ste08l| 
and the references therein, the numerica l analysis of IGA- 
BEM in essentially open. We refer to SBTRI 2 L TMI2 . 
PTG13l| for numerical experiments and to HADI4l| for 
some quadrature analysis. A posteriori error estimation 
has first be en consi dered for Galerkin IGABEM in our re¬ 
cent work FGPI5l| . In the present work, we extend the 
latter result to collocation IGABEM which is preferred in 
practice for its simpler assembly of the stiffness matrix. 
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1.2. Model problem 

Let U C be a Lipschitz domain and F C dfl be 
a compact, piecewise smooth part of the boundary with 
finitely many connected components. For a given right- 
hand side /, we consider the weakly-singular boundary 
integral equation 


Vcl){x) 



y\(l){y)dy = f{x) on F (I.l) 


associated with the 2D Laplacian; see Section [2] below for 
the mathematical setting and the definition of the problem 
related energy norm ||| • |||. With some discrete ansatz space 
Xh C T^(r), the Galerkin BEM computes the unique so¬ 
lution (j)h G Xh of the discrete variational formulation 


Vcjih iphdx 



for all 4ih€ Xh. 


( 1 . 2 ) 


Note that Xh C L^(r) ensures V4>h G C'(r). The colloca¬ 
tion BEM computes (ph & Xh such that 


V(j)h{xj) = f{xj) for all G {xi,... ,XAr^„J, (1.3) 


where the Xj are appropriately chosen collocation points 
with Ncoi := dLim.Xh', see SectionIn either case (I1.2D - 
(ph is computed by solving a linear system of equa¬ 
tions 


1.3. A posteriori error estimation for Galerkin IGABEM 

We assume that Xh is associated to some partition Th 
of F into a set of connected segments. For each vertex z 
of Thj let uJh{z) :=yf{T^Th '■ z & denote the node 
patch. If Xh is sufficiently rich (e.g., Xh contains certain 
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splines or NURBS), it is proved in |FGP15j that Galerkin 
BEM guarantees reliability and efficiency 

Crel^ lll<^ - <Vh-= ^ 

’ (l-4a) 

< Ceff |||<^ - M 

with d:’;j-independent constants CeS^C^ei > 0. Here, Vh := 
f — V4>h denotes the residual and 


vkizr ■■= j [ 

j UJh (z)^ OJ 


\rh{x) - rh{y)\'^ 


uih{z)J bJh{z) 


-y? 


dydx (l-4b) 


is some Sobolev-Slobodeckij seminorm, i.e., the unknown 
BEM energy error is controlled by some computable a pos¬ 
teriori error estimato r rih ■ Estimate (11.41) has first been 
proved by Faermann Fae0nl [ for closed F = dfl and stan¬ 
dard spline spaces Xh based on th e arcleng th parametriza- 
tion. Her result is generalized in FGP15| to a more gen¬ 
eral settin g which also incl udes isogeometric analysis. We 
note that FaeOOl IfGP 15 | show that the efhciency esti¬ 
mate r]h < CeS |||</> — holds even independently of the 
discretization and, in particular, for collocation. 

1.4- A posteriori error estimation for eolloeation IGABEM 

In the present manuscript, we focus on the weighted- 
residual error estimator which has first been proposed in 
CS96I . Car97 | for standard BEM with piecewise polyno¬ 


mials and polygonal F. We prove that for Galerkin IGA¬ 
BEM (11.21) as well as collocation IGABEM (11.31) . there 
holds the upper bound 


C. 


rei 


1 /2 

<dh-=(^J2 d-hiz)^'^ ’ (1.5a) 




1.5. Outline 

Section [2] recalls the functional analytic framework, 
provides the assumptions on F and its parametrization 
7 , introduces the ansatz spaces, and presents an adap¬ 
tive algorithm which is capable to control and adapt the 
multiplicity of the nodes as well as the local mesh-size (Al¬ 
gorithm [221) • Section 13] provides the numerical evidence 
that the proposed adaptive IGABEM is superior to IGA¬ 
BEM with uniform mesh-refinement as well as to adap¬ 
tive standard BEM with piecewise polynomials. Moreover, 
we observe that collocation IGABEM leads to essentially 
the same convergence behavior as Galerkin IGABEM, so 
that an adaptive collocation IGABEM may be favorable 
in pra ctice. Se ction [H recalls the precise statement of (11.41) 
from FGP15l| and gives a proof of (ll.5|) - (ll.7|) . The con¬ 
cluding Section |S] comments on our overall findings, open 
questions, and future research. 

2. Preliminaries 

In this section, we collect the main assumptions on the 
boundary and its discretization and introduce the BEM 
ansatz spaces. Further details on Sobolev spaces and the 
functional analytic setting of weakly-singul ar integral equa ¬ 
tions , are found, e.g., in the monographs HWOSl . IMcLOOI 
SSllj and the references therein. 

Throughout, | • | denotes the absolute value of scalars, 
the Euclidean norm of vectors in , the measure of a set 
in R (e.g., the length of an interval), or the arclength of 
a curve in R^. The respective meaning will be clear from 
the context. 

We write A < B to abbreviate A < cB with some 
constant c > 0 which is clear from the context. Moreover 
A B abbreviates A < B < A. 


with an dA-independent constant Crei > 0. Here, := 
f — V(j)h is again the residual and 

y.h{zf ■■= \^^h{z)\ f \r'^{x)fdx (1.5b) 

J bJhO) 

is a weighted i7^-seminorm, where (•)' denotes the arc- 
length derivative and \ixh{z)\ is the length of the node 
patch. For collocation BEM, we thus control the energy 
error by 

C~s Ih < \H - < C'rel IJ-h (1-6) 

which, however, involves different error estimators. In ad¬ 
dition to the global relation of the error estimators rjh and 
y,h and independently of the discretization, we prove 

Vhiz) < C'ioc/r/t(z) for all vertices z of Th, (1-7) 


2.1. Function spaces 

For any measurable subset w C F resp. any interval 
w C R, denotes the Lebesgue space of all square 

integrable functions with corresponding norm 

■= [ \uix)\^dx. (2.1) 

J OJ 

If M G T^(w) is differentiable along the arc, u' denotes the 
arclength derivative. Define the Sobolev space := 

{u G LA{^) '■ u' G L^(a;)} with corresponding norm 

■= ll^llz, 2 (,^) + (2.2a) 

:= [ \u'{x)\^dx. (2.2b) 

J OJ 

Furthermore, define the Sobolev space := {u G 

LH^) : \\u\\ < oo} with corresponding norm 


where Cioc > 0 depends only on F. 




:= + |n| 



2 

u{yT 


2/1^ 


dy dx. 


(2.3a) 

(2.3b) 
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The dual space of is where duality 

is understood with respect to the extended T^(a;)-scalar 
product, i.e., for u £ and (j) & it holds 

{u\(j))ui= / u{x)(j){x)dx. (2.4) 

J ijj 

We note that C T^(r) C form a Gelfand 

triple and all inclusions are dense and compact. 

Amongst other equivalent definitions of are 

the characterization as trace space of functions in H^(yi) 
as well as equivalent interpolation techniques. All these 
definitions provide the same space but different norms, 
where norm equivalence constants depend only on w. 

2.2. Weakly-singular integral eguation 

The operator V from (HU extends to a linear and con¬ 
tinuous operator V : —>■ with additional 

stability V : L^(r) —>■ iJ^(r). We additionally suppose 
that V is even an elliptic isomorphism, which is satisfied, 
e.g., if diam(fl) < 1. In particular, (I4(-); (•))r is thus a 
scalar product on and the induced energy norm 

uf := {V^P ■ V')r for V' e i?-'/'(r) (2.5) 

is an equivalent norm on 

Given / £ the weakly-singular integral equa¬ 

tion (HU is equivalently reformulated in variational form: 
Find (j) £ i7“^/^(r) such that 

{V(P ; V')r = (/ ; i’)r for all V' e (2.6) 

The Lax-Milgram lemma thus applies and proves existence 
and uniqueness of the solution cp £ i7“^/^(r) of (12.61) 

resp. HU- 

2.3. Parametrization of boundary 

Let r = y. Fi C dV. be decomposed into its finitely 
many connected components F^. Then, 


t £ (a, 6 ] and 7 ^ 0 for t G [a, h). Moreover we assume 

that 7 '^(t) -I- ^ 0 for all c > 0 and t £ [a,b] resp. 

t G (a, b). 

By 7 l : [0, L] —>■ F, we denote the arclength parametri¬ 
zation, i.e., 172(01 = 1 = |72(0Ij its periodic ex¬ 
tension. Then, elementary differential geometry yields bi- 
Lipschitz continuity 


a 


r — 


|7l(s) -7l(0I 


< Cr for all s,t £ 


with 


s — t\ 

^ \s — t\ < J L, for F = on, 
s 7 ^ t G [0,L], for F C Qfi- 


(2.7) 


see, e.g., Ganl4 Lemma 2.1] for the proof for F = dVl 
which even simplifies for F G gj]. Let I C [a, 6]. Suppose 
1^1 < jL for F = i90. Then, (12.71) implies 


Cp O 7L 1771/2(7) < |?7|77i/2(^^(7)) 

< Grju 077,1^1/2(7) 


( 2 . 8 ) 


for all u £ 771/2 (L). 


2 . 4 . Discretization of boundary 

For the discretization, let Th = {Ti,..., T^} be a par¬ 
tition of F into compact and connected segments Tj. The 
endpoints of the elements of Th form the set of nodes 


A/)i 


{zj : j = l,...,n} for F = dil, 

{zj : j = 0,..., n} for F G 


The arclength of each element T £ Th denoted by hr. 
Moreover, the shape regularity constant is defined by 


K(ry :=max{/iT//7T' : T, T'G T^, T n 0}- 


For F = dD, we extend the nodes, elements, and their 
length periodically. Moreover, we suppose 


max hr < |r|/4 for F = dfl. (2-9) 


Il'“ll77i/2(r) — ^2 I 


i) 


see, e.g., [FGP15l Section 2.2]. The usual piecewise poly¬ 
nomial and NURBS basis functions have connected sup¬ 
port and are hence supported by some single F/ each. 
Without loss of generality and to ease the mathematical 
proofs, we may therefore assume that F is connected. All 
results remain valid for non-connected F. 

We assume that either F = dft is parametrized by a 
closed continuous and piecewise two times continuously 
differentiable path 7 : [a, 6 ] —>■ F such that the restriction 
7 l[a,&) is even bijective, or that F G gQ is parametrized by a 
bijective continuous and piecewise two times continuously 
differentiable path 7 : [a, b] —>• F. For F = we denote 
the {b — a)-periodic extension to M also by 7 . For the left 
and right derivative of 7 , we assume that 'y'‘{t) 7 ^ 0 for 


2.5. Discretization of parameter domain 

Given 7 : [a, 6 ] —>■ F, the partition Th induces a par¬ 
tition Th = {Ti,..., TTi} of the parameter domain [a,b]. 
Let a = zo < -Si < • • • < = (/ be the endpoints of the 

elements of Th- We assume T) = [zj-i,Zj], l{Tj) = Tj, 
and '){zj) = Zj. We define 


f^h = 


{zj : j = 1,..., n} for F = dD., 

{zj : j = 0,..., n} for F G 


The length of each f £ Th is denoted by hf. Moreover, 
we define the shape regularity constant on [a, b] by 


K{Th) -maxlhf/hf, : T,T'G 72,7(T’) G 7(T') 7 ^ 0 }. 


Note that K{Th) — K{Th), where the hidden constants de¬ 
pend only on the parametrization 7 . 
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2.6. B-splines and NURBS in the parameter domain 
We consider knots K, := {ti)i^z on R with ti-i < U for 
i £ Z and limi_).±oo ti = ±oo. For the multiplicity of any 
knot ti, we write We denote the corresponding set of 
nodes Jv := {ti : z £ Z} = {zj : j £ Z} with Zj-i < Zj 
for all j £ Z. For i £ Z and p £ No, the z-th B-Spline of 
degree p is defined inductively by 

=3 forp = 0 , 

1 for p > 0 , 


2. 7. NURBS on the boundary 

For F = each node z £ A4 has a multiplicity 
#■2 < P + 1- This induces a sequence of non-decreasing 
knots ICh = (ti)ili on (a, 6 ]. Let Wh = be a se¬ 

quence of weights on these knots. We extend the knot 
sequence (6 — a)-periodically to (ti)igz and the weight se¬ 
quence to (zci)igz by WN+i ■= Wi for i € h. For the 
extended sequences, we also write Kh and Wh- We set 

Nrn]Cu,Wh) :=^^{ICh,yVu)\ia,b)Ol\la]br (2-14) 


where, for t £ R, 


A.p (0 



for ti ti.^p, 
for ti = ti+p. 


We collect some basic properties of B-splines from 

Lemma 2.1 ( dBStil Theorem 6 , Section 2 and page 9-10]). 
For p £ No, the following assertions hold: 


For F C each node z £ A4 has a multiplicity 
<p-\-l such that = # 2 n = p + 1- This induces a 
sequence of non-decreasing knots Kh = {ti)fL_p on [a,h]. 
Let Wh = (wOfci be a sequence of weights. We extend 
the sequences arbitrarily to ICh = (ti)i^z with ti < U+i 
for i G Z, a > ti ^ —oo for z < —p, and b < ti ^ oo for 
i > N, and Wh = {u)i)i^z with Wi > 0. We set 

JpP{ICh,Wh) ■■=^^{JCh,Wh)\iaM°l~"- (2-15) 


(i) Let I = [a, b) be a finite interval. Then, 

{Bi^p\i : z £ Z, i?i^p|/ 7 ^ 0} (2.10) 

is a basis for the space of all right-continuous Jv- 
piecewise polynomials of degree lower or equal p on I 
and which are, at each knot ti, p— ffti times contin¬ 
uously differentiable if p — ffti > 0 . 

(ii) Fori £ Z, Bip vanishes outside the interval [ti_i,ti+p). 
It is positive on the open interval (ti_i,ti+p). 

(iii) For z £ Z, Bi p is completely determined by the p + 2 
knots ti—i ,..., tj^p. 

(iv) The B-splines of degree p form a locally finite parti¬ 
tion of unity, i.e., Bi,p = 1 otz R. □ 


In addition to the knots 1C = {ti)i^z, we consider weights 
W := {wi)ii^z with Wi > 0. For z £ Z and p £ No, we de¬ 
fine the z-th non-uniform rational B-Spline (NURBS) of 
degree p 


r: 


K.W 

'i,p 


■ — ■ — 


'WiBi. 




( 2 . 11 ) 


Note that the denominator is positive and locally finite. 
For any p £ No, we define the vector spaces 


.R'PilC) := 
JFP(t,W) 


aiBfp : Oi £ R 

iGZ 

E ryJC.W 

O-iRi^p ■ Cli G 


. iez 


( 2 . 12 ) 

(2.13) 


Due to Lemma l^TTl dnll-dinll. this definition does not depend 
on how the sequences are extended. 

2.8. Collocation IGABEM 

In this section, we show how to choose the collocation 
points Xj for j = 1 ,..., Ncoi in (IE31). First, we note that 
Lemma [23] © implies that 

{^*,p|[a.b) : Z = 1 - p,..., # 6 -h 1} °l\[a]b) (2-16) 

for F = dll resp. 

{^^*.pl[a,b] : z = 1 - p,..., iV - p} o 7-1 (2.17) 

for F C QGi forms a basis of .yVP(lCh,Wh). Recall ffb = 
p -I- 1 for F C QGl. For simplicity, suppose #6 = p -I- 1 also 
for F = dTl. This gives 


iVeoi = N. (2.18) 

For j = 1,..., the collocation point Xj is defined through 
the arithmetic mean of p -|-2 knots in the parameter domain 

Xj=l(Xj) with Xj := y (2.19) 

p z 


2.9. Adaptive algorithm 


Finally, we recall an adaptive algorithm from our pre¬ 
ceding work FGP15l | , which steers the Iz-refinement of the 
partition 7b as well a s the inc rease of the multiplicity of 
the nodes Nh- While |FGP15| considered ph for Galerkin 
IGABEM, the current focus is on ph and collocation IGA¬ 
BEM. 


Suppose that F is represented by a NURBS curve of 
degree p £ Nq. This induces the initial partition To of F 
with nodes A/q, related nodes A/q in the parameter domain, 
and positive weights Wo- Each node has a multiplicity 
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lower or equal p + 1 , where for F C Qfi or collocation 
IGABEM we suppose #a = #6 = p + 1. For F = dfl, we 
suppose hr < |r|/4 for all T G To- 

As the initial trial space, we consider 

JFp{JCo,Wo) C L^{r) c iT-^/2(F). (2.20) 

Fix an error estimator Qh G {rjinPh}- The nodal con¬ 
tributions Qh{z) from (11.41) resp. (11.51) are used to steer 
knot insertion from ICh to the following knots ICh- The 
new weights W// are uniquely chosen such that the de¬ 
nominator of the NURBS functions does not change. In 
particular, this implies nestedness 


Remark 2.3. (i) While 0 = 1 leads essentially to uniform 
refinement, 0 <C 1 leads to highly adapted partitions. Note 
that the smaller 6, the more iterations of the adaptive loop 
are required. In our experiments below, 6 = 0.75 appeared 
to be an appropriate compromise which led to optimal con¬ 
vergence behavior. 

(ii) The estimate (|2.23l) in step (iv) of the adaptive al¬ 
gorithm can be achieved by some extended ID bisection 
algorithm from \AFF^ lAl . The latter guarantees that the 
overall number of elements is bounded by the sum of ele¬ 
ments in the initial partition plus the number of marked 
elements. □ 


Fi^^{lCh,Wh)QDKP{]CH,WH) 


( 2 . 21 ) 


3. Numerical experiments 


of the related NURBS spaces. Since the weights in Wh 
are just convex combinations of the weights in Wq; it holds 
min Wo < minW ff < maxW// < max Wo. For details, we 
refer to FGP15 . Section 4.2]. 

Then, the adaptive algorithm reads as follows: 


Algorithm 2.2. Input: Adaptivity parameter 0 < 9 < 1, 
polynomial order p G No, initial partition To = 7h with 
knots ICq = ICh, initial weights Wo = Wh- 
Adaptive loop: Iterate the following steps (i)-(vi), until 
Qh is sufficiently small: 

(i) Compute approximation cfh S A^P{lCh, Wh) from Gal- 
erkin BEM (11.21) resp. collocation BEM (11.31) . 

(ii) Compute indicators Qh(z) for all nodes z G Afh. 


In this section, we empirically investigate the performance 
of Algorithm 12.21 for Galerkin as well as collocation IGA- 
BFM in three typical situations: In SectionjXH the bound¬ 
ary F = dH is closed and the solution exhibits a generic 
(i.e., geometry induced) singularity. In Section [331 the so¬ 
lution is smooth on F = dQ, but has certain jumps which 
require discontinuous ansatz functions. In Section 13.41 we 
consider a slit problem. In all examples, the exact solu¬ 
tion is known. This allows to analyze the reliability and 
efficiency of the proposed estimators. 

The boundary part F is parametrised by a NURBS 
curve 7 , i.e., the parametrisation has the special form 

7 {t) = J2c^Rffp'^^t) (3.1) 

iez 


(hi) Determine a set Mh C= Wh of minimal cardinality 
such that 

0gl< Y. ( 2 . 22 ) 

zeMh 

(iv) If both nodes of an element T € Th belong to Aih, T 
will be marked. 


for all t G [a, b\. Here, p G N is the polynomial degree, A-y 
and W^ are knots and weights as in Section l^^ and (C'i)igz 
are control points in which are periodic for closed F = 

dn. 

We choose the same polynomial degree p for our ansatz 
spaces Xh = ^^{iCh, Wh). For the initial knots and weights, 
we choose Kh = Ay, and Wh = Wy. As the ansatz spaces 
are nested, it always holds 


(v) For all other nodes in Mh, the multiplicity will be 
increased if it is smaller than p -\- 1, otherwise the 
elements which contain one of these nodes z € A4h, 
will be marked. 

(vi) Refine all marked elements T £ Th by bisection (in¬ 
sertion of a node with multiplicity one) of the corre¬ 
sponding T £ Th- Use further bisections to guarantee 
that the new partition Th satisfies 

Kifn) < Mfo)- (2.23) 

Update h^ H, i.e., replace Th by Th- 

Output: Adaptively generated partition Th with correspond¬ 
ing solution (fh o,nd error estimator Qh. □ 


l£l2&-W^{th,Wh)\[aM^ (3.2) 

where 71,72 denote the first resp. second component of 
7 . Therefore, this approach reflects the main idea of iso¬ 
geometric analysis, i.e., the same space is used for the ge¬ 
ometry and for the approximation. For adaptive Galerkin 
IGABFM as well as adaptive collocation IGABFM, we 
compare uniform refinement, where AAh = Wh and hence 
all elements are refined, and adaptive rehnement with 9 = 
0.75. In addition, we also consider discontinuous piecewise 
polynomials. Note that this is formally only a special case 
if Wj = 1 for all weights wj of Wh and ffzj = p-\-l for all 
nodes Zj £ Wh ■ 

As basis for the considered ansatz spaces, we use (I2.16P 
resp. (I2.17|) . To calculate the Galerkin matrix, the col¬ 
location matrix, the Faermann error estimator, and the 
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weighted-residual error estimator, we transform the weakly- 
singular integrands into a sum of a smooth part and a log¬ 
arithmically singular part. Then, we use adapted Gauss 
quadrature to comp ute the resulting integrals with appro¬ 
priate accuracy; see Ganl4 Section 5] for details. For the 
weighted-residual error estimator (HID, we replace \uih(z)\ 
by the length | 7 “^(w(z))|, since this eases the calculation. 
Note that \ujh{z)\ ~ | 7 “^(a;?i(z))|, where the hidden con¬ 
stants depend only on the parametrization 7 . 

To calculate, the exact error, we proceed as follows: 
Let (j)^^ € Xh be the Galerkin approximation with 
the corresponding coefficient vector. Let (/)“* € Xh be 
the collocation approximation with the corresponding 
coefficient vector. Let be the Galerkin matrix of 

the /i-th step. With the Galerkin orthogonality and the 
energy norm |||(/)|p = (F^; 0), obtained by Aitken’s A^- 
extrapolation, we can compute the energy error as 




,gal|||2 


= iikiip - ufr 


(3.3) 


resp. 


<(>! 


col III 2 


n 


r 


c‘iip 


(3.4) 


= 11^ - <IIP - {vr\ct - ; icT' - C*))- 


3.1. Laplace-Dirichlet problem 

In the first two examples, we consider the Laplace- 
Dirichlet problem 


—Au = 0 in n and it = gr on T (3.5) 

for given Dirichlet data g € iL^/^(r) and closed bound¬ 
ary r = dfl. The problem is equivalent to the integral 
equation (HU with / = (K + (T)g, i.e. 

V(/)=(K + a-)g on r, (3.6) 

where 

Kg(x) ^ 9 {y)du(y) log(|x - y\) dy (3.7) 

denotes the double-layer integral operator and cr(a:) = 1/2 
for all a; £ r except of the corners, where ij{x) = Q;/(27r) 
with the corresponding interior angle a. The unique so¬ 
lution of (HU is the normal derivative (f) = du/dv of the 
solutio n u G of (13.51) . For more details, see e.g. 

SteOSl . Section 6.3 and 6.6]. 

3.3. Problem with generic singularity 

As first example, we consider the Laplace-Dirichlet prob¬ 
lem (13.51) on the pacman geometry 

O := |r(cos(/3), sin(/3)) ^ 0 ^ ^ | > 


with T = 4/7; see Figure|3Tj The geometry is parametrised 
on [0,1] by a NURBS curve of degree p = 2. We prescribe 
the exact solution of (13.5|) as 

u{x^ y) = r'^ cos (t/3) 


in polar coordinates {x,y) = r(cos/l,sin/3). We consider 
the corresponding integral equation (EH). The normal 
derivative (j) = dujdv of u reads 




( cos(/3) cos (t/ 3) -k sin(/3) sin (t/3)\ , 

\^sin(/3) cos (t/3) — cos(/3) sin (t/ 3)y ’ ^3 


and has a generic singularity at the origin. 



Figure 3.1: Geometry and initial nodes for the experiment from Sec¬ 
tion [SU] 



Figure 3.2: Experiment with singular solution on pacman geome¬ 
try from Section Em The singular solution 0 o 7 is plotted on the 
parameter interval, where 0.5 corresponds to the origin, where 0 is 
singular. 

In Figure EU the solution (j) is plotted over the param¬ 
eter domain. The singularity is located at t = 1/2 and two 
jumps are located at t = 1/3 rep. t = 2/3. 

In Figure 13.31 error and error estimators are plotted. 
All values are plotted in a double logarithmic scale such 
that the experimental convergence rates are visible as the 
slope of the corresponding curves. Since the solution lacks 
regularity, uniform refinement leads to the suboptimal rate 
for the energy error, whereas adaptive refine¬ 
ment leads to the optimal rate In each case. 
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collocation 



Figure 3.3: Experiment with singular solution on pacman geome¬ 
try from Section 13.21 Error and estimator are plotted versus the 
number of knots N. Uniform, ? 7 ^-driven and fi^-dmen refinement is 
considered. 


Figure 3.4: Experiment with singular solution on pacman geometry 
from Section 13.21 The plot shows the efficiency indices |||^^^^||| for 
the estimators G where adaptivity is driven by p^. 


the curves for the two different estimators rjh and jih and 
the error are parallel. In Figure 13.41 we plot the ratios 
rih/l(l>-M resp. ^jih/l(l> - 4>hl- Throughout, these ratios 
stay between 0.5 and 2.7 which underlines an accurate er¬ 
ror estimation for both error estimators. Figure [?75] shows 
the errors of all considered adaptive IGABEM strategies. 
We observe a very similar behaviour. 

For adaptive refinement, Figure llFHl provides a histogram 
of the knots in [a, 6] of the last refinement step for colloca¬ 
tion IGABEM with ph = Ph, for the other adaptive strate¬ 
gies, the output looks similar (not displayed). We see that 
the algorithm mainly refines the mesh around the singular¬ 
ity at t = 1/2. Additionally, the multiplicity at the jump 
points t = 1/3 and t = 2/3 appears to be maximal so that 
the discrete solution (j)h also mimics the discontinuities of 
the exact solution </. 

In Figure 1X71 we finally compare standard BEM with 
discontinuous piecewise polynomials against IGABEM. For 
the error estimation we use the weighted-residual estima¬ 
tor Ph- The output looks similar if rjh is used instead (not 
displayed). All approaches show similar convergence rates, 
however we clearly observe better multiplicative constants 
for Galerkin IGABEM and collocation IGABEM than for 
standard BEM. 



Figure 3.5: Experiment with singular solution on pacman geometry 
from Section 13.21 The errors from all presented adaptive IGABEM 
strategies are plotted versus the number of knots N. 

3.3. Adaptive IGABEM for problem with jump solution 

As second example, we consider the Laplace-Dirichlet 
problem (13.51) on the square = [0,1/2]^; see Figure 
The geometry is parametrised on [0,1] by a NURBS curve 
of degree p = 1. 

We prescribe the exact solution of (Ell) as 
u{x,y) = sinh(27ra;) cos(27ry). 
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Figure 3.6: Experiment with singular solution on pacman geometry 
from Section 13.21 Histogram of number of knots over the parameter 
domain. Knots with maximal multiplicity p + 1 = 3 are marked. 



Figure 3.7: Experiment with singular solution on pacman geometry 
from Section 13.21 The errors from uniform/adaptive BEM with dis¬ 
continuous piecewise polynomials and uniform/adaptive IGABEM 
are plotted versus the number of knots N. 


In Figure [S.lQI we plot error and error estimators. The 
solution ^07 has jumps at the points t = 1/4, t = 1/2, 
t = 3f4 and t = 1 resp. t = 0. As the knots used 
for the parametrisation of T all have multiplicity one, the 
functions of the isogeometric start approximation space 
are continuous at the points t = 1/4, t = 1/2 and t = 3/4. 
Uniform refinement, where only h-refinement takes place, 
leads to the suboptimal rate 0{N~^) for the energy error, 
whereas adaptive refinement increases the knot multiplic¬ 
ity at these problematic points and leads again to the opti¬ 
mal rate 0{N~^^'^). In Figure [IT. Ill we plot the efficiency 
indices ry/j/IIIIII resp. Hh/l\(/-(/hl\- Throughout, these 
ratios stay between 0.1 and 2 . 2 . Figure (?. 12 1 shows the er¬ 
rors of all considered adaptive IGABEM strategies. We 
observe that rjh leads to slightly better results than 
while there appears to be almost no difference between 
Galerkin IGABEM and collocation IGABEM. 

In FigureEHH standard BEM with discontinuous piece- 
wise polynomials is compared against IGABEM. For adap¬ 
tivity, we use the weighted-residual estimator /i/j. The out¬ 
put looks similar if the estimator r]h is used (not displayed). 
We observe that in this example uniform standard BEM 
is superior to uniform IGABEM. This is of course due 
to the fact that standard BEM uses ansatz spaces which 
are discontinuous at the jumps of (p. However, with the 
use of adaptive multiplicity increase this is fixed as can 
be seen in the convergence plot, where we again see that 
adaptive IGABEM leads to better results than adaptive 
standard BEM. It is also interesting that adaptive stan¬ 
dard BEM converges with a better multiplicative constant 
than uniform standard BEM. This is due to the fact that 
the solution is zero on [1/4,1/2] and [3/4,1], wherefore the 
adaptive algorithm uses only few elements in this area. 


7(1) 

7(3/4) 




7(1 

4) 

7(1/2) 



Figure 3.9: Experiment with jump solution on square from Sec- 

Figure 3.8: Geometry and initial nodes for the experiments from tion l3.3l The solution o 7 is plotted on the parameter interval. 
Section 13.31 


We consider the corresponding integral equation dSH). The 
normal derivative p = du/dv of u reads 


<j){x,y) = 2 tt 


f cosh(27ra;) cos(27rj/) 
\^sinh(27rx) cos(27r?/) 


v{x,y). 


It is smooth up to four jumps as can be seen in Figure 


3 . 4 . Adaptive IGABEM for slit problem 

As last example, we consider a crack problem on the slit 
F = [—1,1] X {0}. We parametrize F by a NURBS curve 
of degree p = 1. For f{x, 0) := —x/2 and the single-layer 
operator V, the exact solution of dm reads 


(/)(x,0) 


— X 

\/l — x"^ 





























number of knots N 


collocation 



Figure 3.10: Experiment with jump solution on square from Sec¬ 
tion 13.31 Error and estimator are plotted versus the number of 
knots N. Uniform, 77 ^-driven and p^-driven refinement is consid¬ 
ered. 

Note that 4> S \ L^{T) for all £ > 0 and that cj) has 

singularities at the tips x = ±1. 

In Figure Id.151 error and error estimators for the uni¬ 
form and for the adaptive approach are plotted. The er¬ 
ror is obtained via (1^ resp. (13.41) . where |||(/)|p = 7r/4 
is computed analytically. Since the solution lacks reg¬ 
ularity, uniform refinement leads to the suboptimal rate 
( 9 ( 7 V-i/ 2 )^ whereas adaptive rehnement leads to the op¬ 
timal rate The curves for the two estimators 

and the error are again parallel. In Figure 13.171 we plot 
the efficiency indices 77/^/110III resp. (l)h\l Fig¬ 

ure 13.161 shows the errors of all considered adaptive IGA- 
BEM strategies. Here, jy^i-adaptive Galerkin IGABEM and 
T^h-adaptive collocation IGABEM lead to the best results. 
In Figure 13.181 we compare standard BEM against IGA¬ 
BEM, where we use ph = Ph- While adaptve Galerkin 
IGABEM and adaptive standard BEM lead to optimal 
convergence rates, the best results are achieved with adap¬ 
tive collocation IGABEM. 


Figure 3.11: Experiment with jump solution on square from Sec¬ 
tion [3]3] The plot shows the efficiency indices m for the esti¬ 
mators Ph € where adaptivity is driven by ph- 



Figure 3.12: Experiment with jump solution on square from Sec¬ 
tion [3]3| The errors from all presented adaptive IGABEM strategies 
are plotted versus the number of knots N. 

4. A posteriori error estimation for IGABEM 

4.1. Main results 

For T inductively define the patch a;™(T) C F 

of order m G No by 
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^1{T) :=r, 

W™+1(T) :=\J{r eTh : rnio]^{T)j^(D} 


(4.1) 


















Figure 3.13: Experiment with jump solution on square from Sec- 
tion l3.3| The errors from uniform BEM with discontinuous piecewise 
polynomials and uniform/adaptive IGABEM are plotted versus the 
number of knots N. 


■^O) • • • • • - 9 7(1) 

l^l/S) 7(2/5) 7(3/5) 7(4/5) 


-■- 1 - 1 -'-■- 

-1.5 -1 -0.5 0 0.5 1 1.5 

Figure 3.14: Geometry and initial nodes for the experiment from 
Section 13.41 


The main result of Theorem 14 .2 1 rea uires the following two 
assumptions on Th and Xh for some fixed integer m G Nq: 

(Al) For each T G Th, there exists some fixed function 
ijjT G Xh with connected support supp('07’) such that 

TQsnM^T)Qu:^{T). (4.2) 

(A2) There exists some constant q G (0,1] such that 


collocation 



Figure 3.15: Experiment with singular solution on slit from Sec¬ 
tion niii Error and estimator are plotted versus the number of 
knots N. Uniform, r;^-driven and /i^-driven refinement is consid¬ 
ered. 



Ill - V’T|lt2(,„pp(^^)) < (1 - g) |supp(V't)| (4.3) 
for all T eTh. 


The first theorem shows that these assumptions are, in 
particular, satisfied for NURBS spaces. 


Theorem 4.1 ( [FGPIItI Theorem 4.4]). For p G Nq and 
m := [_p/ 2 ], the space Xh := W/i) satisfies the 

assumptions (A1)-(A2). The constant 0 < g < 1 depends 
only on K(Th), min(>V/t), max(>V;i), p, and 7 . □ 

The main result of |fgpi,^ reads as follows: 


Figure 3.16: Experiment with singular solution on slit from Sec¬ 
tion [Til The errors from all presented adaptive IGABEM strategies 
are plotted versus the number of knots N. 


Theorem 4.2 ( |FGP15 . Theorem 3.1]). For any approx¬ 
imation fih G L^{T), the residual rh = f — Vfih satisfies 
the efficiency estimate 


■= ' < Oeff U 


Vh ■■= 


(4.4) 


zeMh 


with r]h{z) := \rh\H^/^(ojh{z))- If the mesh Th and the dis- 
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collocation 




number of knots N 


The following two theorems are the mathematical con¬ 
tributions of this work to the field of IGABEM. They apply 
to both, Galerkin IGABEM as well as collocation IGA¬ 
BEM. 

Theorem 4.3. For any approximation (ph £ 

Th ■= f - V(j)h, the indicator r]h{z) := is 

bounded above by the weighted-residual indicator fJ-hiz) ■= 

Phiz) <V2Cr ti-hiz), (4.6) 

where Cr > 0 is the constant from ism . 

If collocation IGABEM as in Section 12.81 is used, the 
patch ujP^^{T) contains a collocation point and therefore 
a root of the residual r^, for each T G Th- Hence, the 
condition of the following theorem is fulfilled with m = 
p-\-l. 

Theorem 4.4. Suppose that wither (fh £ is the Galerkin 
solution of (EH), where Xh satisfies (A1)-(A2), or that the 
residual ru = f — Vcfh has at least one root in each w™(r) 
for all T gTh (ind some fixed m £ No- Then, 

1 /2 

E ^^hiz)) (4.7) 

z^Afh 


Figure 3.17: Experiment with singular solution on slit from Sec¬ 
tion [3]4] The plot shows the efficiency indices m for the esti¬ 
mators ph G where adaptivity is driven by p^. 



withpLhiz) := \u!hiz)\'^^‘^ IWhWL'^iuihiz))- The constant Cre\ > 
0 depends only on E, m, k(7/i), and, in the first case, q. 


4-2. Proof of Theorem^ 

We o nly need the following lemma, whose proof is 
inspired by NPVllL Proposition 2.2], where an analogous 
assertion for norms instead of s eminor ms is found. The 
assertion itself is also stated in |GE01I . Lemma 7.4] in a 
more general way. Indeed a similar version of (14.81) hold s 
even for the i7®-seminorm, 0 < s < 1. However, in [CFOlj, 
the proof is only given for the hardest case 1/2 < s < 1. 


Lemma 4.5. For any connected w C P, whose length sat¬ 
isfies jw] < jL j/P = dil, there holds 


|w|||u'||i 2 (^) for alluG H\r). (4.8) 


Figure 3.18: Experiment with singular solution on slit from Sec¬ 
tion [33| The errors from uniform BEM with discontinuous piecewise 
polynomials and uniform/adaptive IGABEM are plotted versus the 
number of knots N. 


Proof. We recall that for a finite interval / C M, H^{I) 
coincides with the space of all absolutel y contin uous func¬ 
tions on / with derivative; see, e.g., Evald , page 306]. 
Step 1: First we consider I = (0,1) and prove 


Crete space Xh satisfy assumptions (A1)-(A2), the Galerkin 
solution (fh £ Xh of EH also satisfies the reliability esti¬ 
mate 


klifi/2(7) < (4-9) 

We use the transformation theorem, with r = pis — t) +t 
and s — t = tr, as well as the Cauchy Schwarz inequality 


III/' - < Grel?7/»- (4-5) 


The constant CeS > 0 depends only on P, while C^ei > 0 
additionally depends on m, ^(7/), and q. □ 
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to get 


lifi/2(7) 


J j u[s)-u(t) 

a 


ds dt 


s — t 

f(0,s) - /(O.t) ^'(^) 


s — t 


ds dt 


hh +^) 


ds dt 


< 


IJI JI 


u'i^p{s — t) + t)\ dpdsdt 

' / / \u'[pa + dp da dt. 

I J(-t,i-t) Ji 


We formally extend u' by zero to R. This and the Fubini 
theorem lead to 


1^1771/2(7) — 
< 


I / \u'{pa+ t)\^ dpda dt 

JIJ (.-1,1) JI 

/ / / \u'[pa + t)\^ dt da dp 

Ji J(.-i,i) 7r 

hi 11 = 


'(- 1 , 1 ) 

Step 2: If / C R is an arbitrary finite interval, it holds 

1111^1/2(7) ^ 2|/||u|^i(j). (4-10) 

Without loss of generality, let I = [c,d) be open. We 
define the function //(gp) : (0,1) —>■ R : r !->■ u[T[d — c) + 
c). Obviously, it holds //(op) S i?i(0,l) with ^^(t) = 
[d — c)u' [T[d — c) + c). The transformation theorem with 
s = a[d — c) + c, t = T[d — c) + c, and r = p[d — c) + c, 
and (14.91) yield 


1771/2(7) 


ll 


u[s) — u[t) 


s — t 


ds dt 


f f 

u[a[d — c) -I- c) — u[T[d — c) + c) 

J( 0 , 1 ) Ji 0 ,l) 

a — T 


ds dt 


— |ii(o,i)li/i/ 2 ( 04 ) < 2|m(o,i) 


= 2 


[ \u[op)[p)\^dp = 2 [d-c) [ \u'[r)\'^dr 

'(0,1) Ji 


= 2|/||u|/i-i(/). 

Step 3: We show (14.81) . Let / be a real interval with 
7 z,(/) = w. Then, (|2.8D and (|4.10l) give 


|ll|77i/2(a;) = ^ 1l\h^/‘2(I) 

< 2Cp|a;||M o = 2 Cp|w|||u'||p 2 (,^). 

This concludes the proof. □ 

4-3. Proof of Theorem 

We use t he following estimate from FaeOOl . Lemma 2.3]; 


see 


Ganl4 Proposition 2.13] for a detailed proof. 


Lemma 4.6. There exists a constant Ci > 0 such that, 
for all u € it holds 

I|ll|l77i/2(r) - X] |l''l77i/2(a;;,(z)) + C*! X/ I|ll|li 2 (r) • 

zeW TeTfe 

(4.11) 

The eonstant only depends on F and K,[Th)- D 

Proof of Theorem EZl If the residual is orthogonal to some 
Xh satisfying (A1)-(A2), the assertion follows at once from 
Theorem 14.21 in combination with Equation (14.6F If the 
residual has local roots, we first note that 


\\4’ - 4’h\\H-i/2(Y) — 11/ - VcjlhWm/^iT) — l|l’/i||77i/2(r)i 

(4.12) 


since V is an isomorphism. The hidden constants only 
depend on F. 

Taking u = rh in Lemma 14.61 it only remains to es¬ 
timate the sum Wl^t)- Note that shape 

regularity yields \ujff[T)\ < (2m-I-1 )k(7)i)'"/it- Replacing 
T by uj'^[T), we apply Friedrich’s inequality to see 


/l,p^||r;^||p2(p) < ^T^II^'^llL2(tj™(T)) 

T&Th TGTh 


TGTh 


K[Tr 

hx 


h{urn^T)) 


TGTh 

<[2m + lfi,[Thf^ Y, hT\K\\l^^^r) 


TGTh 


<{2m + ifi.[rHf^ Y KWIIk;.lli2(.(z))- 

zGXh 


This concludes the proof. 


□ 


5. Conclusion 


5.1. Analytical results 

In this work, we considered adaptive BEM for weakly- 
singular integral equations Vcf = f associated to elliptic 
PDEs in 2D. As model example served the 2D Laplacian, 
but the results apply as long as V : H~^^'^[T) -T iL^/^(r) is 
an elliptic isomorphism. With the residual rh '.= f — Vcfh, 
we transferred the weighted-residual error estimator 


Ph = ll/i^/V;,l|p2(p) 


(5.1) 


proposed in CS96l . [^ar97l] from standard BEM with lowest- 
order polynomials to IGABEM, where we considered the 
Galerkin method as well as collocation. Eor either dis¬ 
cretization, we proved that ph is reliable 


11^ - < GrelM^; (5.2) 
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see Theorem 14.41 In our preceding work FGP15| , we con¬ 
sidered the residual error estimator 


r]h = 


(E / [ 

^ \r. Ju)h.{z)Juj 


\rh{x) - 


dydx'j (5-3) 




-y? 


proposed in |Faennl| . In |FGP1,^ , we transferred this esti¬ 
mator from standard BEM with piecewise polynomials to 
IGABEM. Independently of the discretization, we proved 
the general efficiency estimate 


yh^CefilcI)- (5.4) 

while our proof of the converse estimate ||| ()) — ||| < Crei Vh 
is restricted to Galerkin IGABEM. However, the combina¬ 
tion of (15.21) and (15.41) provides also full error control 

C~^ T]h < \14>- 4>hl\ < Crel y-h (5.5) 

for collocation IGABEM computations in 2D. Moreover, 
this estimate implies the global relation rjh y-h-, and we 
even proved 

yh{z) < C\oc rih{z) for all z G Afh (5.6) 

for the respective nodal contributions defined in (El resp. 
El; see Theorem 14.31 which holds independently of the 
discretization employed. 


5.2. Numerical results 

We proposed an adaptive algorithm which is capable to 
steer the mesh-refinement as well as the knot multiplicity 
in Galerkin and collocation IGABEM computations; see 
Algorithm l2.21 Numerical experiments in Section [3] under¬ 
line that generic singularities of the (unknown) exact solu¬ 
tions lead to reduced experimental convergence behavior if 
the underlying mesh is not appropriately graded. This is 
a well-known fact for standard BEM with piecewise poly¬ 
nomials, but also applies to IGABEM. Consequently, the 
gain of adaptive IGABEM (resp. the loss in case of uni¬ 
form meshes) is huge due to the higher-order ansatz func¬ 
tions of IGABEM, and therefore adaptivity seems to be a 
must to exploit the full potential of isogeometric analysis. 
In several numerical experiments, we showed that the pro¬ 
posed algorithm is capable to recover the optimal order of 
convergence. The gain of IGABEM is that the algorithm 
chooses smooth NURBS, where the exact solution appears 
to be smooth, while discontinuities and singularities are 
well detected and appropriately resolved. Compared to 
standard BEM with discontinuous piecewise polynomials, 
this leads to a smaller number of degrees of freedom for 
comparable accuracies. 

For collocation IGABEM as well as Galerkin IGABEM 
and independently of the (uniform or adaptive) mesh-re¬ 
finement, we observed that 

Vh ^ \U - Ml ^ yh, (5.7) 


i.e., both error estimators are efficient and reliable. The 
efficiency indices r]h/\l4> - Ml and yh/ll4> - Ml appeared 
to be < 3, i.e., the overestimation of the energy error is 
very moderate. We note that only the equivalence ijh — 
|||</> — for Galerkin IGABEM as well as the bounds 
lll</> - Ml ^ yh and r]h < |||<(> - M have thoroughly been 
proved mathematically. 


5.3. Open questions and future work 

All considered numerical experiments show optimal con¬ 
vergence of the estimator and the error. Unde rstanding 
this observation mathematically in the spirit of CFPP14 j 
is one of our goals for future research. However, it is ques¬ 
tionable if an analogous version of the re duction prop¬ 
erty on refined element domains GFPPIdl (A2)] can be 
proved for the Faermann estimator rjh. Indeed, this is 
yet an open problem even for sta ndard BEM with piece- 
wise polynomials; see FFME+I3| . where at least conver¬ 
gence of an /i-adaptive algorithm with ph is analyzed. For 
the weigh ted-residual error estimator yh the axioms of 


CFPP14j are satisfied for standard Galerkin BEM with 
piecewise polynomials, see CFPP141 Section 5.4]. For 
collocation IGABEM there remain two challenging math¬ 
ematical questions: First, one needs further investigation 
on the unique solvabilit y of the d iscrete system. Second, 
the quasi-orthogonality CFPPl4 (A3)] is unclear for col¬ 
location methods. 

As mentioned, we observed in all numerical experi¬ 
ments reliability as well as efficiency of the used error es¬ 
timators. However, it remains to mathematically verify 
the reliability estimate ||| </) — <(>;i ||| ^ r]h for collocation BEM 
and the efficiency estimate yu ^ |||<(> — fihf + osc, at least 
for some higher-order oscillation terms osc. Again, these 
estimates are yet open problems even for standard BEM. 
For lowest -order Ga lerkin BEM, the efficiency estimate is 
proved in AFF+131 . Theorem 4] under additional regular¬ 
ity assumptions on the Dirichlet data g in (13.5|) . 

Finally, the ultimate goal is of course to analyze and 
apply the estimators ph and yt in 3D Galerki n IGABE M. 
For 3D one has to c onsider e.g., T-splines SSE+I3l| or 
hierarchical B-splines BGI5l| , because, in contrast to mul¬ 
tivariate N URBS, they naturally allow for local mesh re¬ 
finement. Fae02l| shows that gu is reliable and efficient 


for stan dard BEM with piecewise polynomi als, whereas 
GMSOl j proves reliability for yh. In |gFPPi 4 Section 5.4] 


optimal convergence of adaptive /i-refinement for yh is 
proved. The estimate ph fy yh as well as plain conver - 
gence for 77 /i-based adaptivity is analyzed in [FFMEjJHj. 
The transfer of the mentioned results from standard BEM 
to adaptive IGABEM leaves interesting and challenging 
questions for future research. 
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